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ABSTRACT 

Results of an extensive investigation of probability distribution functions 
(pdfs) for Rayleigh-Benard convection, in the hard turbulence regime, is pre- 
sented. It is seen that the pdfs exhibit a high degree of internal universality. 
In certain cases this universality is established within two Kolmogorov scales 
of a boundary. A discussion of the factors leading to universality is presented. 
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1 Introduction 


Probability distribution functions (pdfs) long have been of importance in 
discussing turbulent flows (see Pope 1 for a comprehensive review to 1985). 
They serve as a valuable tool in presenting data and give valuable insights 
into range of behavior. Recent work 2-6 has placed a new emphasis on their 
study in turbulent flows. Direct numerical simulations 6-9 and improved ex- 
perimental techniques 10-12 have permitted more detailed pdfs, and allow us 
to go well beyond such low order moments as variances < 7 , skewness s and flat- 
ness /, common in the statistical description of turbulent flows. Part of the 
recent interest in pdfs is due to the experiments of the Chicago group, 10-12 
and in particular to the exponential tails they observed for the pdfs of tem- 
perature fluctuations. The wide tails of such pdfs underline the essential role 
that intermittency plays in turbulence phenomena. However, earlier investi- 
gations also led to such exponential distributions 13-15 . 

On the theoretical side, Sinai and Yakhot 2 considered the pdf of a passive 
scalar and made remarkable progress in describing the limiting distribution 
function. Subsequent to this, Yakhot and coworkers 3 ’ 16 extended this ap- 
proach and presented arguments leading to exponential pdfs for temperature 
fluctuations and vorticity. Kraichnan has also made a notable extension to 
this approach. 5 More recently Kraichnan 17,18 has produced a relatively sim- 
ple closure model based on the mechanisms at work in fluid flow which lead 
to wide skirted pdfs. 

In this paper we present the results of an extensive study of pdfs as 
generated in a computational investigation of the Rayleigh-B6nard (R-B) 
convection problem. The present results are based on a substantially longer 
simulation than was considered in our earlier publications. 6,7 For purposes of 
comparison we estimate the turnover time by 

( 1 ) t = H/qrm, 
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where H is the height of the computational cell and q rmi represents the rms 
speed. In such terms our earlier work was based on 4.65 turnover times while 
the present study is based on 41.7 turnover times. 

The specification of the problem as well as the details of the calculation 
appear in our earlier work 6,7 . In brief, we simulate turbulent thermal con- 
vection as described by the Boussinesq equations 19 . The vertical coordinate 
is denoted by z and the horizontal by x and y. The horizontal planform is a 
square and the aspect ratio (width to height) is 2\/2. The Rayleigh number 
is .98 x 10 4 times the critical value, which places the simulation in the range 
of hard turbulence in the terminology of the Chicago group. 10-12 The velocity 
satisfies slip boundary conditions on the impermeable horizontal boundaries 
and temperature is specified on these walls. Periodic boundary conditions 
are imposed in both horizontal directions. 

The calculation uses 96 equal grid spacing in the vertical direction, with 
a like number in each of the horizontal directions. Since the Nusselt number, 
Nu, is roughly 23 for this calculation and since 

< 2 > 

where 6 Is the thickness of the thermal sub-layer, two grid spacings lie in the 
sublayer. (The sub-layer thickness is also a measure of the Komogorov scale.) 
These few facts are useful in interpreting the material that follows since 
various quantities will be specified by their vertical gird spacing locations. 

As will be seen the pdfs exhibit a remarkable degree of universality. By 
universality we mean the unanticipated lack of dependence on strongly in- 
homogeneous directions. As discussed in Section 6 the same universality 
is expected to extend to other convection problems if the pdfs are suitably 
normalized. 
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2 Preliminary Considerations 

In general if <x represents a dependent variable, then we will denote by 



the pdf which will be plotted in our figures. The bar signifies that the quan- 
tity has been averaged in time as well as over a horizontal plane. As previ- 
ously reported for Benard convection 6 ’ 7 temperature pdfs vary considerably 
in the vertical direction, Figure 1 shows the pdf of temperature fluctuation, 
T', at six different elevations. (See the cited reference for more details of 
the calculation and notation.) The temperature pdf at the edge of the ther- 
mal sub-layer ( z = 2/96) though not symmetric, has skewness close to zero. 
Above the bottom sub-layer (z = 7/96,16/96, and 26/96) the pdfs exhibit 
a peak corresponding to a small negative value and a long positive tail, re- 
sulting in significant positive skewness. This indicates the existence of hot 
thermal plumes near the bottom boundary and at the mid-plane the pdf is 
symmetric. Figure 2 shows the comparable pdfs for the horizontal velocity 
u and vert iced velocity w. The parabolic shape of the u pdf in this log-linear 
plot indicates the Gaussian nature of the distribution. At the mid-plane the 
pdf of w exhibits a similar Gaussian nature, whereas below the midplane the 
pdf of w is positively skewed similar to the temperature pdfs but to a much 
lesser extent. 

In our earlier papers we confirmed the existence of an exponential range 
to the pdf for fluctuations in temperature in the midplane as first exhibited 
by the Chicago group. This result was based on less than three log units 
range in the pdf. With the increased database, we now have reliable data for 
a range of almost six log units as exhibited in Figure 1. While the pdf in the 
midplane is still fit by an exponential over the first two to three log units the 
overall picture is quite different. The reason for the sharp down turn in the 
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tails of the pdf is easy to understand. Since the temperature fluctuations are 
bounded below and above by the wall temperatures, the skirts of the pdf, 
P(T'), must terminate at finite values. This tendency is clearly exhibited 
by the midplane pdf shown in Figure 1. A brief summary describing the 
shape of the P(T') at the midplane is: (1) first, analytical considerations 
dictate that the pdf must be rounded at the center; (2) next, this is followed 
by wide skirts, fit by an exponential, and which indicate relatively large 
fluctuations; (3) finally, limited temperature fluctuations dictated by the 
wall temperatures force a rapid falloff of the skirts. 

In viewing the preceding pdfs and those to be presented later it is impor- 
tant to keep in min the symmetries of the problem. This is especially true 
of the pdfs measured in the midplane of the cell z = 48/96, since some pdfs 
become symmetric only in the midplane. As has been shown, the problem 
under discussion has a sixteen-fold symmetry group 20 . We have not made 
use of this to extend the database (and thus sharpen the pdf curves) but in- 
stead have verified these symmetries to support the correctness of data. (E.g. 
P(u), not shown, is virtually identical to P(u), Figure 2(a).) In considering 
the pdfs it is useful to split these into two classes. Those that must have 
symmetric distributions as a result of an inherent symmetry and those for 
which there is no k priori symmetry requirement. In general, if we consider 
a pdf P(d), and if there exists an admissible transformation of the group, G, 
such that 


(4) Ga = — a 
then this implies that 

(5) P(&) = P(—a). 

Consider for example the horizontal velocity u. Under refection in the 
plane x = 0,u — * — u, and hence P(u) is symmetric for all z values. This 
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is clearly indicated in Figure 2a. Alternately, both the vertical velocity w, 
and the temperature fluctuation T' go into their negatives under reflection 
in the midplane, but not elsewhere. Thus their mid-plane pdfs should be 
symmetric and as seen in Figures 1 and 2(b), this is the case. However no 
group operation produces (4) at other z elevations and there is no reason 
for these pdfs to be symmetric for other values of z. This is borne out by 
Figures 1 and 2(b). 

We can carry this line of reasoning somewhat further and consider w t and 
T x , neither of which have a group operation leading to (4). Figures 3a and 3b 
show the pdfs for each of these quantities. As can be seen the corresponding 
pdfs show no symmetry at any z values. This can be given a physical inter- 
pretation. Away from the two thin sub-layers near the top and the bottom 
boundaries the pdfs are negatively skewed. To see this this consider a parcel 
of fluid located at some height z' travelling up. Since the vertical velocity 
and temperature fluctuations are well correlated, the probability that this 
parcel of fluid has positive temperature fluctuation is high. After a small 
time interval the fluid parcel will be at a higher elevation z 1 + 8z' surrounded 
by relatively colder fluid. Therefore the probability that this parcel will ac- 
celerate up ( w z > 0) and contribute larger temperature fluctuation (T' z > 0) 
is high. Similarly if we consider a parcel of fluid with negative temperature 
fluctuation moving down, the probability that the parcel will accelerate down 
(w x > 0) with larger temperature fluctuation (T z > 0) is large. Therefore 
both w z and T t show a peak in their probability distribution at a (small) 
positive value and this peak is compensated by a less steeper negative tail, 
in order to yield zero mean values. For the P(w z ) this negative tail can be 
associated with the less probable event of rapid deceleration of both the up 
moving cold parcel of fluid and the down moving warm parcel of fluid. In 
other words, mild acceleration up or down is more probable than mild de- 
celeration, but rapid deceleration (up or down) is more probable than rapid 
acceleration. Similar interpretations can be given for the gradient of the tem- 
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perature fluctuations as well, and in general the pdf of a quantity provides 
us with a view of its range of behavior. 

For comparison with the pdfs of Figure 3, we consider U t which has a sym- 
metry leading to (5) for all z. The corresponding pdf is shown in Figure 4. 
This clearly indicates the symmetry but in addition shows a remarkable de- 
gree of universality in the ^-direction. Such universality will be encountered 
for the majority of the pdfs which will be presented below. Further discussion 
of this unexpected property will be presented in the following sections. 
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3 Probability Distribution Functions 


As already seen in Figure 4 the pdf of u z (or equivalently v x ) exhibits a 
surprising universality across the convective cell. To explore this property 
further we present, in Figure 5, the components of vorticity. The pdfs for the 
horizontal components flj and 0 2 are virtually identical, only Qi is shown, 
and in addition exhibit a high degree of universality. Vertical vorticity, 0 3 , on 
the other hand shows a significant departure from universality. All vorticity 
components lead to pdfs which satisfy (4), either through reflection in a plane 
of constant x or y. 

Given that mechanisms are at work that force universality, we should 
expect in general that boundaries will produce significant departures from 
universality. To contrast the results of the pdfs for the horizontal and vertical 
components we recall that 



It follows from the boundary conditions for this particular R-B convection 
problem 19 that = 0 (and 0 2 = 0) at a wall, while 0 3 is unrestricted at a 
wall. As a general principle one should expect that all permissible fluctuations 
will appear in the corresponding pdf with an appropriate probability. This 
can be seen in Figure 5. Thus which is unrestricted at the boundary has 
a wide skirt at z = 2/96, while which is restricted at the wall can only 
develop a slight skirt in the neighborhood of the wall z = 2/96. As a result 
universality does not appear in P(fl 3 ) until we approach the central region, 
where the pdf is fit by the universal curve. On the other hand universality is 
quickly achieved for P(0 1 ) and P(0 2 ). It should be noted that in the central 
region the pdf of all three components of vorticity are well fit by the same 
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exponential over at least three decades. Further the pdfs of those velocity 
derivatives that constitute the three components of vorticity (u v , u z ,v x ,v z ,w x 
and w y ), coincide with the above exponential distribution (for example see 
Figure 4). 

We have already seen in Figures 1 and 3, examples of pdfs that show 
a strong departure from universality. For each of these pdfs the quantity 
in question does not satisfy the symmetry property (4) except possibly at 
the center plane. Both u x and v y are linked to w z through the continuity 
equation. Since neither of these quantities have a transformation leading to 
(5), their pdfs lack universality. The pdfs for u x and v y closely resemble 
those shown in Figure 3(a) and are not presented here. (Only the pdf of w z 
at z = 2/96 is significantly different.) [Although u x and v y are identically 
distributed, they are correlated, since the convolution of P(u x ) with itself 
does not produce P(u/ Z ).] 

The pdfs for w x , with high accuracy lie on those of Oj , Figure 5(a). By 
contrast with the pdf in w, Figure 2(b), this is universal except very near 
the wall. Thus taking a derivative of w, which then leads to a symmetrizing 
transformation, (4), produces a quick transition to universality. In addition, 
the value of w x is pinned to be equal to zero at both the top and bottom 
boundaries. This transition to universality is less true for T x , the pdf of which 
is shown in Figure 6. We recall, however, that the pdf of T' itself, Figure 1, 
is far less universal and significantly more skewed than the pdf of w , Figure 
2(b). 

In general it can be observed that temperature statistics obey universality 
over a narrower region near the midplane than their velocity counterpart. 
This reluctancy towards universality can be, at least partially, attributed to 
the highly intermittent nature of temperature and its derivative signals. For 
example, at the midplane while the pdf of u and w are Gaussian the pdf of 
T' is exponential and while the pdf of vorticity is exponential the pdf of T' x 
is flatter than an exponential. 
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4 PDFS for Higher Derivatives 

By taking derivatives of flow variables we emphasize the smaller scales. 
It is a widely held view that small scales forget their large scales origins. 
This concept lies at the heart of the -5/3 Kolmogorov range 21 and of the 
exponential dissipative range 22 . Both these energy ranges represent univer- 
sal behavior of turbulent flows. Thus, as derivatives of flow quantities are 
increased their pdfs depend more heavily on the higher wavenumbers and a 
tendency towards universality should be expected. Pdfs of derivatives say 
something about distributions in small eddies, the higher the derivative the 
smaller the dominant eddy. What must be regarded as remarkable is that 
even pdfs of first derivatives show universality. It would appear that more 
than the above mentioned universal ranges are at work in producing univer- 
sality, and very probably that dynamics is important in establishing at least 
some of the universal features. 

In this section we consider the effect of taking additional derivatives of 
flow quantities on the corresponding pdfs. In Figures 7-9 we exhibit a selec- 
tion of second derivatives of flow quantities. These show a higher degree of 
universality than their first derivative counterparts. Before commenting on 
these in detail, we observe that all the pdfs shown in these figures exhibit a 
flared out skirt and are no longer fit by a simple exponential. To account for 
this effect, we refer to an argument given sometime ago by Kraichnan 23 which 
demonstrates that greater intermittency is to be expected as the wavenum- 
ber is increased. His argument is largely independent of Reynolds number. 
Thus in viewing pdfs of higher derivatives of flow quantities we should expect 
increasing intermittency and hence wider flaring skirts in the pdfs. 

In Figures 7(a) and 7(b) we show the pdfs for u zz and u yt both of which 
exhibit strong universality. On the other hand pdfs for u yy shown in Figure 
7(c) exhibit some departure from universality. u yt is identically equal to zero 
at both boundaries, due to the stress-free boundary condition, therefore a 
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strong tendency towards universality can be expected. Though the bound- 
ary condition for the other three quantities are not prescribed, variations in 
u zz will be small since it is the vertical gradient of u x which is identically 
equal to zero at boundaries. As seen earlier, functions that are fixed at the 
boundaries quickly attain universality and within the universal regime the 
vertical gradients as well as the horizontal gradients can be expected to be 
bounded and exhibit strong universality. On the contrary, both and u xy 
are neither fixed at the boundary nor do they represent the derivative of a 
quantity which exhibit universality. 

We do not show the three pdfs corresponding to w xx ,w zz and w XZ} since 
they virtually lie on the universal curve of Figure 7(b). The pdfs of the first 
two have no transformation under which (5) holds. Symmetry is nevertheless 
established. This confirms that skewness is primarily a property of large scale 
structures. 9 Both w xx and w zz ( by continuity) are restricted to be zero at the 
boundary, whereas w xz is the vertical derivative of the universal function w x . 
Therefore all three pdfs show universality and to excellent approximation fall 
on the same curve, Figure 7(b). In Figure 8 we show pdfs for T xx and T xz . 
The pdfs for the last quantity does not have a symmetrizing transformation 
leading to (4). Nevertheless except near the wall for P(T XZ ) the pdfs are 
symmetric and exhibit a tendency towards universality. This is even more 
true for T XI which also does not have a symmetrizing transformation but has 
pdfs very well fit by Figure 8(a). 

Thus we conclude that whatever the mechanism forcing universality, as 
might be expected, it is more effective on higher derivatives. 


10 



5 Other Mean Quantities 


Figure 9 contains plots of T^^Urm^Wrm,, as well as < WT 1 > (and 
its covariance) and gj < T' >. In certain instances for a range of z, the 
curves are well fit by power laws and this has been indicated. Arguments 
leading to power law dependence are similar to those leading to the inertial 
sublayer (log-layer) for turbulent wall bounded flows. 24 In brief once away 
from the diffusive layers adjacent to a wall the only available length scale 
is the distance to the boundary. Prandtl 25 applied this to the convection 
problem and one finds 24,26 


( 8 ) 


T rm. « * 


-1/3 


W r 


oc z 


1/3. 


dz 


oc z 


-4/3 


in this region. T^, and w rm , are the rms temperature and vertical velocity 
fluctuations (Kraichnan, 27 also using mixing length theory, has considered in 
detail the case in which the Pr number is allowed a wide range of values.) 
A simple argument which leads to (8) is that since in the inviscid region the 
convective heat transport 


(9) = H 0 

is a constant, this implies that w and T' scale reciprocally with z. Thus, in 
the vertical momentum equation the lead terms are 

dw i 

( 10 ) w a7 ccT 

from which the scaling in (8) follows. 

While T rmt does show a sensible range for a (-l/3)-power, the same is 
not true for w rmi . Though the convective heat transport < wT' > and (the 
covariance < wT 1 > /u; rmj Tj m4 ) are constant across the layer (.20 < z < .35) 
except near the boundary, w does not scale reciprocally since the region 
over which T' oc z 1 ^ 3 (.05 < 2 < .2) falls closer to the boundaries. We also 
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note that dT' fdz } does not follow a -4/3-power. Townsend 28 ' 29 observed 
a -2 power, indicated as a dashed line in the Figure, a value also given 
by Carroll 30 and in a computation by Eidson et al. 31 This value was first 
proposed by Malkus 32 from theoretical arguments. However we point out 
that a -9/4 power, the continuous straight line of the Figure, is a better fit 
to the calculation. 

The lack of a consistent universal scaling regime is reminiscent of the 
situation for boundary layer flows, where an argument similar to that lead- 
ing to (8) produces the log layer for the mean velocity 33 . However as has 
been known for some time, not all quantities follow the universal scaling. 
To account for the lack of universality Townsend 34 and later Bradshaw 35 
postulated the notion of active and inactive portions of the flow. In their 
explanation they invoke the idea of large (integral scale) eddies entering the 
boundary regions at infrequent times and marring the universality of some 
quantities. 

In the present case, the underlying scaling arguments assume a semi- 
infinite domain in the vertical direction without horizontal wind, 26 and so 
might be thought to be applicable outside the dissipative sublayer at a wall. 
Recent experiments by the Chicago group indicate a sustained symmetry 
breaking wind. Also for the numerical simulations already discussed 6,7 , the 
roll motions act as a winds since their time scale is of relatively long duration 
although they have zero mean. In any event this more general case changes 
the above dimensional reasoning and is treated in Monin and Yaglom 26 who 
find general classes of possible functional dependences. We do not pursue 
this further since scaling does not appear to be a dominant effect. 

Since all the pdfs considered are given in the normalization, (3), it is of 
interest to present the remaining variances of the quantities considered with 
respect to their variation in the vertical direction, z. The rms fluctuation of 
the three components of vorticity are plotted in Figure 10. Figures 11, 12 
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and 13 show the rms fluctuation for the first and second derivative quantities. 
Since pdfs of the raw variables follow from these plots this completes the 
single point probabilistic description of all the relevant quantities. 

While scaling does not appear to account for the universality another 
feature of the pdfs does help explain the collapse of pdfs onto universal 
curves. We observe that the above symmetric pdfs are well approximated by 
the family of curves 

(11) P{a) = Ce~^ 1 


where C = C(p, k ) is easily chosen so that there is unit area under the curve. 
p = 2 corresponds to a Gaussian distribution, p = 1 corresponds to an 
exponential distribution and p less than one indicates a flatter intermittent 
distribution. All the even moments of this distribution depend on p and k, 
while the odd moments are identically equal to zero. Closer to the boundaries 
the distributions are intermittent and the pdfs are flared out corresponding 
to p less than 1. As universality is approached the value of p decreases 
and in the universal regime since the level of intermittency can be expected 
to be uniform in the value of p i.e. it will be independent of the vertical 
location. Although k, which measures the width or the standard deviation 
of the pdf, can vary in the universal regime for the raw variables, once the 
raw variables are normalized by the rms fluctuation the value of k is dictated 


to be 




TP / 2 


L^UlJ 


so that the normalized standard deviation is unity. Therefore 
in the universal regime p and k (which now depends only on p) are constants. 
Therefore as long as the level of intermittency (or p) reaches an asymptotic 
value, indicating universality, the collapse of the normalized pdfs can be 
explained. As already mentioned the approach to universality is rapid in the 
case of variables that are fixed at the boundaries. 
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6 Concluding Remarks 


An overall survey of the pdfs that have been exhibited leads one to believe 
that there is an active mechanism forcing the pdfs toward universal form. 
Both the degree to which this is true and the resulting shape of the pdf 
depend on the dominant relevant eddy size i.e. on the number of derivatives 
being considered. In general the smaller the eddies, the more quickly is 
universality established, and the more intermittent (flared skirts) the shape 
of the pdf. As has been observed by She et al 9 skewness is essentially a large 
scale property, while flatness is a small scale effect. 

As might be expected the presence of boundaries mars universality and 
wide departures from universal behavior can be expected in the neighborhood 
of a boundary. A boundary can act in two extremely different ways. It can 
pin down a fluctuation. E.g. w, T' ,w X) u z , . . . are all restricted to vanish at a 
boundary. On the other hand fi 3 , T' x , w z , u v , u x , . . . are all unrestricted in that 
any fluctuation in those quantities is permitted at the boundary. Another 
important ingredient in determining the form of a pdf is symmetry. If a 
transformation of the form (4) is applicable then a symmetric pdf results. In 
certain instances this is only achieved at the midplane, e.g. w and T' . While 
in other cases symmetric pdfs must be obtained (assuming that there exists 
sufficient data) at all elevations, e.g. u, w z ,u z , — 

As a general rule, at a boundary, and in its neighborhood, if there is no 
restriction imposed on the quantity, one should expect all manner of possible 
fluctuations to appear. This helps to explain why for example, P(Q 3 ) has 
relatively wide skirts in the neighborhood of the boundary, and as a result 
shows a slower tendency to universality. By comparison u z must be zero 
at a boundary and thus does not develop wide skirts in the neighborhood 
of a boundary. Unlike P((l 3 ),P(u z ) tends very quickly to a universal form, 
Figure 4. 

The cases of P(T') and P(w) are of interest to consider from the perspec- 
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tive of the present discussion. In both instances the quantity in question is 
forced to vanish at the boundary. This however has little effect in restricting 
the corresponding pdf. Although each must be symmetric in the midplane, 
this is not true elsewhere and the passage to the symmetric pdf is difficult 
to characterize. 

The mechanism which is responsible for the tendency towards universal- 
ity is not obvious. Kraichnan 5 ’ 18 has produced a simple heuristic model of 
intermittency based on a closure approximation which exhibits independence 
of Re. This may be indicative of the processes at work and which bring about 
universality. Further investigation of this effect is clearly indicated. 

While our deliberations are based on the computation of R-B convection, 
at one value of Ra, it seems clear that they should be generalize to other 
flow geometries and other values of the control parameter. The approach to 
universality which we have followed in terms of the vertical distance from 
a wall is doubtless better expressed in terms of a wall normalized variable. 
For channel or boundary layer flows this would be the usual wall normal 
coordinate, while in the present case z/5 is the suitable variable. When the 
pdfs are expressed in such terms we anticipate that universality will also 
hold with varying control parameters. A study of available channel flow data 
strongly indicates universal forms for the corresponding pdfs. 37 . 
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Figure 1: Probability distribution function of normalized temperature fluc- 
tuation T' = r '~ T ' at six different heights from the bottom boundary. 

V(T'-T') w# . . 

z = 2/96 is at the edge of the bottom thermal sub-layer, z = 7/96 is m the 
region of hot thermal plumes and z — 16/96 is roughly at the edge of the 
plume region, z = 26/96,37/96 and 48/96 are in the turbulent core. The 
results corresponding to the top half can be extended from the bottom half 
based on symmetry conditions. 
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Figure 3a: Same as Fig. 1 but for w~ t . 
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Figure 4: Same as Fig. 1 but for 
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Figure 5a: Same as Fig. 1 but for Hi. 
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Figure 5b: Same as Fig. 1 but for CI 3 . 
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Figure 7a: Same as Fig. 1 but for u«. 
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Figure 7c: Same as Fig. 1 but for Uyy. 
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Figure 8b: Same as Fig. 1 but for T zt . 
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Figure 11: Rms fluctuation in u x ,u z ,w x and w z against vertical distance. 
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Figure 12: Rms fluctuation in T X) T' Z ,T' XX and T' zz against vertical distance. 
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Figure 13: Rms fluctuation in u w ,u tz and u yz against vertical distance. 
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